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Propagating left/right asymmetry in the zebrafish embryo: one-dimensional 
model 
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ABSTRACT During embryonic development in vertebrates, left-right (L/R) asymmetry is reliably generated by a conserved 
mechanism: a L/R asymmetric signal is transmitted from the embryonic node to other parts of the embryo by the L/R asymmetric 
expression and diffusion of the TGF-/3 related proteins Nodal and Lefty via propagating gene expression fronts in the lateral 
plate mesoderm (LPM) and midline. In zebrafish embryos, Nodal and Lefty expression can only occur along 3 narrow stripes 
that express the co-receptor one-eyed pinhead (oep): Nodal along stripes in the left and right LPM, and Lefty along the 
midline. In wild-type embryos, Nodal is only expressed in the left LPM but not the right, because of inhibition by Lefty from the 
midline; however, bilateral Nodal expression occurs in loss-of-handedness mutants. A two-dimensional model of the zebrafish 
embryo predicts this loss of L/R asymmetry in oep mutants [15]. In this paper, we simplify this two-dimensional picture to a 
one-dimensional model of Nodal and Lefty front propagation along the oep-expressing stripes. We represent Nodal and Lefty 
production by step functions that turn on when a linear function of Nodal and Lefty densities crosses a threshold. We do a 
parameter exploration of front propagation behavior, and find the existence of pinned intervals, along which the linear function 
underlying production is pinned to the threshold. Finally, we find parameter regimes for which spatially uniform oscillating 
solutions are possible. 
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1. Introduction 

In vertebrates, organs such as the heart and brain develop in 
an invariant left-right (L/R) asymmetric fashion, with situs 
inversus, or complete organ mirror-reversal, being a rare oc- 
currence [1, 2, 3], The best-understood mechanism by which 
L/R bias is first established is that special cilia around the 
embyronic node drive a leftwards fluid flow [4]. Disrupting, 
or reversing this flow results in randomization [5], or reversal 
[6] of organ asymmetry. The cilia appear to be conserved in 
vertebrates [7], but it is controversial whether the nodal flow 
mechanism is decisive in all vertebrates, since asymmetry 
of some other origin (involving electrical potentials due to 
polarized gap junctions) is observable before the cilia even 
begin moving, and this may determine organ asymmetry in 
Xenopus [8], 

Our concern in this paper is with the subsequent propa- 
gation of this L/R signal from the node to the rest of the 
embryo. It is first transmitted to the lateral plate mesoderm 
(LPM), resulting in L/R asymmetric expression of the TGF- 
(3 related signaling molecules Nodal and Lefty [9, 10]. Nodal 
is expressed exclusively in the left LPM, and Lefty in the 
midline, along gene expression fronts that propagate from 
the posterior to the anterior of the embryo from approxi- 
mately the 12 to 20 somite stages [11]. This propagation de- 
pends on the diffusion of Nodal and Lefty, as well as their 
production rates which are promoted or inhibited by these 



same molecules. Subsequently, downstream genes (such as 
Pitx2) are activated and direct L/R asymmetric organogene- 
sis [12]. 

Nodal and Lefty constitute a classic activator-inhibitor 
system that is highly conserved among vertebrates [1]. Such 
a model was introduced and studied for the initiation of 
Nodal and Lefty expression in the LPM and midline, respec- 
tively, at the level of the node in the mouse [13]. This was 
spatialy one-dimensional, representing the L/R axis, and did 
not attempt to model the longitudinal propagation of the L/R 
signal. 

In this paper, we extend the ideas in this model to the ze- 
brafish embryo, in which Nodal and Lefty gene expression 
only occurs in cells carrying the co-receptor one-eyed pin- 
head (oep) found on narrow stripes 3-5 cells wide in the 
midline and LPM [14]. (Note that the Nodal analogue in the 
zebrafish is called "Southpaw", but we will call it "Nodal" 
for consistency with the terminology of Ref. [13]). An ad- 
vantage of this system is that it is quasi-one-dimensional in 
the longitudinal direction, because the most important vari- 
ables are Nodal and Lefty production on these 3 stripes, and 
the propagation of their gene expression fronts is clearly uni- 
directional in the posterior-to-anterior direction. 
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A two-dimensional model of the Nodal/Lefty system in 
the zebrafish embryo has been set up and simulated [14, 
15], as we review in Section 2. In Sections 3 and 4, we 
work through the consequences of a one-dimensional ideal- 
ization of it, representing not the L/R axis (as in [13]) but the 
anterior-posterior axis. We idealize Nodal and Lefty produc- 
tion to have step-function turn-ons, when a threshold func- 
tion linear in their concentrations crosses a threshold. We 
classify all possible behaviors in the limit that Nodal is un- 
affected by Lefty, and demonstrate the existence of inter- 
vals along which the Lefty profile must be pinned to that 
of Nodal. In Section 5, we find parameter regimes that allow 
uniform spatial oscillations to occur. 



2. Two-dimensional model 

In this section, we outline what would be a minimal model 
for the two-dimensional Nodal/Lefty system. This will mo- 
tivate the simplification we make, in section 3, to a one di- 
mensional model of similar form. The degrees of freedom in 
the model we use here are only the concentrations N(r, t) of 
Nodal and L(r, t) of Lefty, as found unbound in the extracel- 
lular fluid. In this section we let r represent a two-component 
position (x, y) on the embryo, which is practically a flat two- 
dimensional layer enveloping the yolk cell; in later sections 
on the one-dimensional model it will be replaced by the vari- 
able y. Also, t is the time. 

The geometry is idealized so that the three oep-expressing 
stripes extend indefinitely for y > 0, and are straight, paral- 
lel, and of unvarying widths, as shown in Fig. 1. That means 
that, away from the baseline at y = 0, the model system has a 
translational invariance. That permits one to talk in a mathe- 
matically precise way about a limiting behavior in the model, 
in particular uniform motion of a front. Whether this is per- 
tinent to the actual embryo depends on whether initial tran- 
sients persist for a short time compared to the stages in which 
the L/R signal propagates, and for a short distance compared 
to the inter-stripe spacing or to the embryo's length, which 
we know to be the case [11]. 

Lefty is produced at a rate s i (r, t) [concentration per unit 
time per unit area] which can be nonzero only at points ly- 
ing within the midline stripe, since only those cells can pro- 
duce Lefty. The function Si(r, t) depends, in our simplified 
model, only on the local, instantaneous concentrations of 
Nodal and Lefty; the exact functional dependence is spec- 
ified in Sec. 2.2. More precisely, the production function de- 
pending on N(r, t) and L(r, t) should rather be the transcrip- 
tion rate of the Lefty mRNA, while the source rate sl (r, t) 
of Lefty protein is in turn proportional to the mRNA concen- 
tration, and thus should lag the production function we use 
by roughly the mRNA lifetime. (Indeed, a model including 
that lag was successfully fitted to quantitative data in [14], 
Chapter 3.) But in this paper, we assume that the source rate 
is proportional to the instantaneous production function; that 
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FIGURE 1 Layout of the 2D problem [15]. Nodal is only ex- 
pressed in the dark shaded regions, and Lefty is only expressed 
in the light shaded region. Initially the left tailbud expresses 
more Nodal than the right. 



is equivalent to assuming the Lefty mRNA is transcribed co- 
piously and decays quickly compared to other time scales. 

Similarly, Nodal is produced at a rate sn {r, t) which is 
nonzero only within stripes in the left LPM and right LPM, 
being described by a production function of the same form 
as that of Lefty in the midline, and with the same production 
function in the right and left LPM, so no prior asymmetry is 
assumed. As with Lefty, we have made an approximation of 
short-lived mRNA so that the mRNA content is not treated 
as an independent variable. 

In the case of Nodal, there is a second place where SN(r,t) 
is nonzero, which are the two tailbuds, located as shown in 
Figure 1. This production is not modulated by N(r, t) or by 
L(r,t), and may be assumed constant in time for simplic- 
ity (in reality it gradually diminishes [14]). The tailbud pro- 
duction has some left-right asymmetry, which in our model 
is responsible for all downstream asymmetries. The origin 
of the asymmetry is imagined to be in the nearby Kupffer's 
vesicle [16] (the organ in which the special cilia are found in 
zebrafish). 



2.1 Diffusion behavior 

Once produced, Nodal and Lefty are assumed to be passive, 
noninteracting densities that satisfy a diffusion equation with 
diffusion constants and Dl respectively. Furthermore, 



they are assumed to have degradation rates, respectively r 



N 
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and . Although one expects the two signaling molecules 
to have similar values of these parameters, we allow them to 
be different in general, because (as we derive in the results 
sections) this inequality induces qualitiatively different be- 
haviors. (The differential equations incorporating these pa- 
rameters are exactly like Eqs. (3.2a) and (3.2b) of the one- 
dimensional model, except the second derivative with re- 
spect to y is replaced by the two-dimensional divergence.) 

Characteristic length and velocity scales can be constructed 
from the diffusion constant and degradation times: 



In=VD 



nt~n; 



TN 



(2.1) 



(2.2) 



and we define II and vl similarly. These scales determine 
the typical range in space that the signal molecules travel, 
and the typical speed of the production fronts we shall study. 



2.2 Conditions for Nodal and Lefty production 

A key property determining the behavior is that expression 
(and production) of either Nodal or Lefty is promoted by the 
presence of Nodal and inhibited by the presence of Lefty. 
Thus, both production functions sjv (r, t) and sl (r, t) are in- 
creasing as a function of N and decreasing as a function of 
L. As we explain next, one can consider several versions of 
our model with different functional forms for the production 
functions. 

Let us first lay out how production is actually regulated. 
Along the LPM and midline stripes, Nodal and Lefty pro- 
tein bind to the plasma membrane with the help of oep co- 
receptors [14]; we may say "oep" is in one of three states: 
unbound, Nodal-bound or Lefty -bound. We assume that only 
a negligible fraction of all Nodal and Lefty produced are thus 
sequestered by the receptors, so that uptake or release of the 
signaling molecules is not included in Sjv(r, t) and Si(r, t). 

In turn, Nodal-bound oep induces a signaling cascade that 
promotes transcription of Nodal (in the LPM) or Lefty (on 
the midline). The fraction of Nodal-bound oep, <f>(N, L), is 



(KN,L) 



N 



N + L 



(2.3) 



The oep receptors on the midline are believed to be identical 
to those on the LPM stripes, so we use the same function 
4>(N, L) in both stripes with the same equilibrium constant 
k (where fc _1 has the units of concentration). 

The Nodal production function s n is assumed to be pro- 
portional to a Hill function 



i 



i + {(t>* N /<f>) nH 



(2.4) 



This is a rounded step from to 1 centered at a threshold pa- 
rameter c/)*at; the step gets sharper as the Hill exponent tih 



grows large. A microscopic rationalization for the form in 
Eq. (2.4) would be that the Nodal gene's regulatory region 
cooperatively binds njj copies of the signal or transcription 
factor induced by Nodal-bound oep. (Ref. [13] used a differ- 
ent production function than (2.4) that also approximates a 
step function with a parameter controlling the sharpness.) 
The actual production rate is 



SN(r,t) = s N0 f N ((j>(N,L)) 



(2.5) 



where snq is the maximum, or saturated, production rate. A 
similar function Jl (<fi) , with (it is expected) a rather different 
threshold 4>*l, and a saturated production rate slq, controls 
Lefty production. 

Part of the motivation for including the intermediate func- 
tion <f> in the model is to represent mutants in which oep is 
under- or over-expressed. That would have the effect of mul- 
tiplying cj) in (2.4) by a constant, or equivalently of dividing 
</>*jv and <j>*L by that constant [14]. 

We adopt two simplifying assumptions, either of which 
could be made independently of the other. The first pertains 
to the threshold for turning production on or off. In place 
of the nonlinear <j> function (2.3), we follow Ref. [13] and 
instead use two threshold functions linear in the concentra- 
tions: 

C N (N,L) = C NN N - C NL L (2.6a) 



C L (N,L) = C LN N-C LL L 



(2.6b) 



where N = N(r,t) and L = L(r,t) are Nodal and Lefty 
concentration, respectively, and Cnn, Cnl, Cln, and Cll 
are constants. Production of N or L respectively turns on 
when this function exceeds a threshold parameter Cn* or 
Cx*, which should be positive, since a certain concentration 
of Nodal is needed to turn on either production, even in the 
absence of Lefty. We set Cn* = Cl* = 1 without loss of 
generality by scaling the four coefficients in Eqs. (2.6). 

The relation of Eqs. (2.6) to the <f> function in (2.3) is that 
the condition to be over threshold, <j> > (f>*N or <j> > <p*L, is 
(by inserting Eq. (2.3)) mathematically equivalent to 



l)N-kL > 1; 



*(£ri - 1)W - kL > 1. 



(2.7a) 
(2.7b) 



Clearly, we can (roughly) identify the coefficients in (2.7) 
with those in (2.6), in particular this viewpoint implies Cn l = 
Cll- 

Note that although the condition to be at threshold is cap- 
tured exactly by the linear functions, away from threshold 
the production will depend on N and L in a different way 
when Tin < oo. The linear form of Eqs. (2.6a) and (2.6b) 
will be more convenient mathematically. 

The second simplification we make is to assume njj — > 
oo in (2.4), i.e. a step-function turn-on of the production. In 
this case, Nodal production s n (r, t) and Lefty production 
si,(r, t) are given by: 
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and 



SN(r,t) — , C N < 1 

s$ s N (r,t) < s N0 , C N = 1 

SN(r,t) = s N0 , C N > 1; 

'« L (r,t) = , C£<1 

0<« £ (r,t) <« i0 , C £ = l 

s L (r,t) = s L0 , C L > 1. 



(2.8a) 



(2.8b) 



2.3 Behavior of the two-dimensional model and 
comparison with experiments 

The model initializes with zero Nodal or Lefty concentra- 
tion in the domain (Fig. 1), and with the 2 tailbuds produc- 
ing Nodal that diffuses to the base of the 3 stripes. The left 
tailbud produces more than the right, which is responsible 
for all downstream asymmetries. Successful initiation at the 
base of the stripes requires that the Nodal produced in the 
tailbuds is sufficient to start self-sustaining Nodal produc- 
tion in the left LPM, but not the right LPM; production be- 
gins when thresholds given in (2.8) are exceeded. The times 
at which this occurs differs for all 3 stripes depending on the 
difference between 0*jv and the amounts of Nodal pro- 
duced by the left and right tailbuds, as well as the position of 
the tailbuds (closer to or farther from the midline, in a L/R 
symmetric way). Lefty, once produced on the midline, dif- 
fuses to the left and right LPM in equal amounts. Thus, the 
threshold function (2.6a) is constrained such that Nodal pro- 
duction turns off in the right LPM due to Lefty inhibition, but 
stays on in the left because of greater N there. This occurs 
in simulations for widely-varying parameter sets, for differ- 
ences in tailbud production rates of 30%. As this asymmetry 
is decreased, the coefficients in (2.6a), as well as the geome- 
try of the problem, must be tuned more finely for successful 
initiation to occur. 

Experiments in zebrafish find that after an initial phase the 
wavefronts of Nodal and Lefty gene expression proceed at a 
fixed speed [14, 11], and indeed this is reproduced in simula- 
tions: once initiated in a self-sustaining way, Nodal produc- 
tion on the left LPM switches on from the base upwards, and 
the leading edge of production moves at an asymptotically 
fixed rate along the stripe, followed by a front of Lefty pro- 
duction, moving at the same asymptotic rate on the midline. 

The main role of Nodal and Lefty in this model is to trans- 
mit the L/R asymmetric signal from the tailbud. Lefty in 
particular acts as a "barrier" preventing initiation of Nodal 
production in the right stripe: the Lefty concentration there 
is sufficiently high that Nodal diffusing from the left LPM 
is insufficient to initiate production there. Another neces- 
sary condition for successful propagation is that the function 
(2.6b) must not fall below 1 along the midline, or else Lefty 
production will cease, and Nodal produced in the left LPM 
will not be prevented from diffusing to the right LPM and 



initiating production further up the right stripe, leading to 
loss of L/R asymmetric propagation. Indeed, in Lefty knock- 
out mutants, bilateral Nodal expression is observed [11, 13]. 

We intended to model the space and time-dependent be- 
haviors in this model to, ultimately, 

(i) confirm the correctness of the basic picture by its qualita- 
tive agreement with experiments; 

(ii) explain, and predict, the phenotypes of various mutations 
involving this system; 

(iii) characterize the robustness of how the initial bias in the 
tailbud is amplified (as is manifested in the error rate); 

(iv) show how the numerical values of the various parame- 
ters can be inferred (or as least bounded) on the basis of 
experiments. 

However, while it is understood in principle how we get a 
moving front, the mathematical formulation involves convo- 
lutions in two dimensions and does not include closed forms 
of simple functions. Furthermore, simulations suggested var- 
ious possible regimes in which production might turn off 
again after the initial passage of a pulse, or in which the 
midline Lefty production had a gradual onset in space, much 
less sharp than the Nodal front. In order to analytically rep- 
resent the functional form, and to comprehensively classify 
all regimes of the asymptotic behavior, we turned to a one 
dimensional model, which is the main focus of the rest of 
this paper. 

3. One-dimensional model 

The 2D model may be simplified to ID by making the fol- 
lowing simplifications. Since Nodal and Lefty expression 
occur on narrow, parallel stripes, we may replace the finite 
width of each stripe by a ID line. Also, the concentration 
profile on one stripe is transmitted (with a diffusive lag) to 
another through the intervening space; this space can be re- 
placed by lines as well, and we end up with a 5-line ID 
model with diffusion between stripes given by 

Ni (y, t) =k D N i _ 1 (y,t- At) + k D N t+1 (y, t - At) + ... 

(3.1) 

where N{ and iVt±i are Nodal concentrations on adjacent 
stripes, and kfj and At describe the rate of diffusion. 

In this paper, we will only focus on the front propagation 
of Nodal and Lefty gene expression in the steady state, in 
which case we further neglect inter-stripe diffusion and thus 
simplify the system to ID. 



3.1 Diffusion equations 
The partial differential equations for Nodal and Lefty are: 



8N(y,t) n d 2 N(y,t) 1 . . , . 

' =D N — N(y,t) + s N (y,t) (3.2a) 



at 



dy 1 



TN 
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dL(y,t) d 2 L(y,t) 1 

= d l — — L(y,t) 



Of 



dy 1 



TL 



s L (y,t) (3.2b) 



The equations appear linear at first glance, but the s jy and 
sl terms are in fact nonlinear, since (according to (2.8a) and 
(2.8b)) they depend on (N, L) via step functions. Neverthe- 
less, if either Cn{v) < 1 or CV(y) > 1 throughout an in- 
terval, then sjv(y) is constant (either or sno) in that inter- 
val, in which case the differential equation for N(y) may be 
solved, and similarly for L(y). Evidently, the key parameter 
in such a solution is the position y at which the threshold is 
crossed and production turns on or off. 

Surprisingly, it is also generically possible that the thresh- 
old function is pinned at 1 throughout a finite interval (see 
Sec. 4.2). Within such an interval, the production function 
is indefinite (the middle case in Eqs. (2.8a) and (2.8b)) and 
the diffusion equation can no longer be used to find the solu- 
tion. However, the threshold condition becomes an equality 
interval and gives the solution in such intervals. The compli- 
cations in our results have to do mostly with handling such 
"pinned" intervals. 

The simplest possible solution is for constant, homoge- 
neous production, e.g. sn (y , t) = sno for all (y, t) for Nodal. 
The solution is N(y, t) = N , or L(y, t) = L in the analo- 
gous case for Lefty, where 



= TnSnq', 



L — t l slo- 



(3.3a) 



(3.3b) 



The middle of any interval of production looks locally like a 
piece of this homogeneous system, so it is not surprising that 
Nq and Lq serve as a reference level for all solutions. In par- 
ticular, since sn(u, t) ^ sno everywhere, any solution must 
have N(y, t) ^ Nq everywhere, and similarly for L(y, t). 



3.2 The traveling wavefront solution 

Experiments show the Nodal front is typically ahead of the 
Lefty front [11, 14]. If it is sufficiently far ahead, the Lefty 
concentration there is negligible, and it is sufficient to con- 
sider a Nodal-only situation, in which the front of Nodal pro- 
duction moves at constant speed v: 



SN(y,t) 



s NQ , y <vt 
, y > vt. 



(3.4) 



The front shows a constant shape to a viewer traveling at 
velocity v, i.e. N(y,t) = N(y — vt). When this is inserted 
into Eq. (3.2a), it becomes an ordinary, second-order, linear 
differential equation, which is solved by an exponential: 

(N 0[ l- 9N e^] ,y-vt<0 

{N (l - g N )e~ KN{y - vt) ,y-vt>0. 

(3.5) 



Substituting this ansatz into Eq. (3.2a), we find that kn and 



are given by 



"TV 



V 


1 


L>n 


' N 


v , 


1 




" w 



= 



(3.6a) 



(3.6b) 



By matching the boundary conditions at y — vt = 0, we find 
g N is given by 

9n = -^-r (3-7) 



K N 



Essentially, this means that the shape of the concentra- 
tion profile depends on the relative magnitudes of v and ujy- 
Physically, for v > and v 3> vn, the wavefront is trav- 
eling forward too fast for there to be much diffusion from 
the producing region to the non-producing region, and so the 
concentration at the front is low (g^ — > 1); conversely, for 
v < and \v\ ^> v^, the concentration at the front is high 
for the same reason. (It is reminiscent of the Doppler effect, 
except that the Nodal signal spreads by the diffusion equa- 
tion rather than the wave equation.) 

There is in fact a symmetry relating solutions with v > 
and v < 0, whereby the profiles of an advancing front of 
speed v and a retreating front (reflected in y) of speed —v add 
up to a uniform profile, because adding the production of an 
advancing front and that of a retreating front simply obtains a 
uniform producing line. Furthermore, when v = the height 
of the front is simply ^Nq, and the stationary concentration 
profile simplifies to: 



\N e-i»/ lN 



,y<o 
,y>o 



(3.8) 



(Notice that we set y = to be the position of the front.) 
Since the coefficients in the threshold equations (2.6a) and 
(2.6b) are in general different, the Lefty wavefront is in gen- 
eral displaced by a distance Ay from the Nodal one. We have 



L(y-vt) 



' L [l - g L e KL( - v - vt+Ay) ] 
Ln(l - g L )e- K '^ v - vt+Ay) 



y — vt < —Ay 
y — vt > —Ay 
(3.9) 

where Ay > signifies that the Nodal front is ahead of the 
Lefty front, and kl, k' l and g^ are defined analogously to 
Eqs. (3.6) and (3.7). Finally, to find v and Ay, we observe 
that N = N (l - 5 at) at the Nodal front, and L = Lq(1 - 
gh) at the Lefty front, and substitute (3.5) and (3.9) into the 
threshold functions (2.6); assuming that Ay > 0: 

C NN N (l - g N ) - C NL L (1 - g L )e-< Ay = 1 (3.10a) 

C LN N [1 - g N e- K » Ay ] - C LL L (1 - g L ) = 1 (3.10b) 
from which we can solve (non-trivially) for v and Ay. 
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3.3 Nondimensionalization 

Certain combinations of parameter changes are trivial, in the 
sense that the solutions look the same apart from rescalings 
of the distance, time, or concentrations. Since we are faced 
with the difficulty of exploring a large parameter space, we 
wish to discover the minimum number of nontrivial param- 
eters. To this end, we will scale all variables and parameters 
in the problem so as to make them dimensionless. 

First, we introduce the scaled Nodal and Lefty concentra- 
tions 

N(y) = ^, (3.11a) 



No' 
L 

^0 



(3.11b) 



so that N(y), L(y) «C 1. 

Correspondingly, Cjviv. C^l, C l n and C LL are simply 
the coefficients in the (scaled form of the) threshold func- 
tions (2.6a) and (2.6b): 



Cnn = CnnNo; 
Cnl = CnlLo'- 
Cln = ClnNq; 
Cll = CllLq- 



(3.12a) 
(3.12b) 
(3.12c) 
(3.1 2d) 



For example, Nodal production turns on when C^^N — 
CnlL > 1, and similarly for L. The parameters (3.12) are 
pertinent even in a spatially uniform situation. 

The remaining parameters relate to time or length scales. 
We scale length and time such that the parameters for Nodal 
become unity: 



II = Il/In; 
vl = v L /v N ; 
v = v/vn- 



(3.13a) 
(3.13b) 
(3.13c) 



The last two parameters are relevant (in a moving steady 
state) if and only if v 0. Implicitly, Eqs. (3.13) also nondi- 
mensionalize the time scale: 



ft = Il/vl 



(3.14) 



In total, we have seven nontrivial parameters for the one- 
dimensional problem, defined in Eqs. (3.12) and (3.13). By 
rescaling space and time with y — > ^- and t — > we find 
that Eqs. (3.2) become 

dN(y,t) _ d 2 N(y,t) 
dt dy 2 

d~L{y,t) - . d 2 L(y,t) 



- N(y,t) + s N (y,t) (3.15a) 



Ilvl- 



—L(y,t) + s L (y,t) (3.15b) 



dt dy 2 
where ^ sjv(y, t) ^ 1 and < §L(y, t) ^ — 
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FIGURE 2 Possible behaviors on the (JV, L) phase plane. 
Green and red lines represent N(N , L) = and L(N, L) = 0, re- 
spectively. Stable (unstable) fixed points are marked by closed 
(open) circles; the circle in (d) with a dot is a fixed point that may 
either be stable or may become unstable to a cycle, as elabo- 
rated in Sec. 5. The topology is classified into cases (a),(b),(c), 
and (d), following Sec. 3.4, according to whether the N isocline 
is to the right or the left of the L isocline, or how they cross (if 
they do). Each case has several subcases according to where 
the respective isoclines interect the bounding square. Cases (a) 
and (c) are compatible with Eqs. (2.7), which requires that both 
isoclines have the same vertical intercept; cases (a),(b), or (c) 
might be compatible with the condition Cnl = adopted in 
Sec. 4, which required that the N isocline is vertical. Either of 
cases (a) or (b) is consistent with a fixed point at < Z* < 1, 
associated with "pinned" intervals throughout which (N, L) is 
exactly on theX = isocline (Sec. 4.2). 

3.4 The (iV, L) plane: zero dimensional case 

In Sec. 3.2, asymptotically in either direction both TV and 
L approach constant, uniform values. Therefore, in order 
to classify possible wavefront solutions, it is helpful first 
to classify all possible uniform solutions. Fig. 2 shows the 
(TV, L) phase plane with lines representing the N — and 

L = isoclines. The most important feature of these dia- 
grams is the relation of the two isoclines to each other and to 
the bounding lines: this determines the possible steady states. 

There are 4 qualitative cases for the relation of N to L: 

1 . iV-isocline is to left of L-isocline 

2. L-isocline crosses iV-isocline left to right 

3. i-isocline is to left of iV-isocline 

4. iV-isocline crosses L-isocline left to right. 

(There are other mathematical cases in which one of the iso- 
clines does not pass through the square at all, but that is ob- 



Biophysical Journal: Biophysical Letters 



L06 



Biophysical Journal: Biophysical Letters 



viously to be discarded, since there would be no way ever to 
produce the corresponding signal.) 

These phase planes are generalizations of the wild-type 
scenario: as noted after Eq. (2.7), if we justify the linear 
threshold equations (2.6) from the single-binding site recep- 
tor behavior (2.3), we necessarily find Cnl = Cll', the ge- 
ometrical expresion of this on the (N, L) plane is that the 
two isoclines cross on the negative L axis, hence only cases 
(a) and (c) could be realized. However, our purpose here is to 
study the general behavior of these equations that could arise 
in any embryonic system with propagating Nodal and Lefty 
fronts (since this is conserved in all vertebrates), and perhaps 
even mathematically equivalent equations having a quite dif- 
ferent biological interpretation. In general, binding of Nodal 
and Lefty might be cooperative at the receptor; furthermore, 
gene regulation further downstream might be cooperative. 
Thus, the general form of the isoclines is nonlinear and we 
expect that all these topologies are imaginable. 

There is always a stable fixed point at N = L = 0; this 
is the only fixed point in case (iii), which is thus trivial. If 
the iV-isocline and the L-isocline both intersect the upper 
border, that is if Cnn — Cnl > 1 and Cln — Cll > 1, we 
have a stable fixed point with both N and L saturated, which 
underlies the basic wavefront behavior of Sec. 3.2 above. 

However, if the L-isocline intersects the right edge, to the 
right of the TV-isocline, i.e. when Cln — Cll < 1 as may 
happen in either case (i) or case (ii), then the stable fixed 
point at N = 1 exists at a less than saturated value, which 
we will define i„: 



1. ee ^L_l ( 3.i 6 ) 

Cll 



It is not obvious that L being pinned at L* < 1 exists in 
the case with spatial variation, but we will in fact show in 
Sec. 4.2 that analogous "pinned" intervals arise in ID travel- 
ing wave solutions. 

Finally, in case (iv), there is a fixed point at which both N 
and L are less than saturated; this may be stable, but also may 
be unstable to periodic oscillations around the fixed point, as 
elaborated in Sec. 5. 

4. Results: steady-state fronts 

To explore the large parameter space of the ID model, we 
will consider special cases; the aim is a taxonomy of the pos- 
sible behaviors of the model, with the hope that any general 
case will be qualitatively similar to one of these behaviors 
seen in the special cases. In particular, we concern ourselves 
with the steady state limiting behaviors, after all transients 
have died out. Any nonzero Nodal and Lefty wavefronts will 
then be traveling at identical speed v. 
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4.1 C NL = and v = 

For the largest part of our story, we consider the case that 
Cnl = 0, that is, Nodal is unaffected by Lefty (i.e. Nodal is 
an autonomous variable). We first solve for the behaviors of 
the one-component system in which Nodal activates Nodal, 
which already has moving front solutions; then, we treat the 
Nodal concentration N(y, t) as if it were externally imposed 
and solve for another one-component problem representing 
the Lefty concentration field. The reasons for choosing what 
seems to be a major simplification are (i) the reverse approx- 
imation, in which Lefty is autonomous, has only the trivial 
solution with Lefty off (since Lefty only inhibits); (ii) exper- 
imentally, the Nodal front leads and the Lefty front follows; 
thus, it is plausible that the inhibition from Lefty has no im- 
portant effect on the Nodal front (serving only to prevent the 
initiation of a front on the other side); (iii) in Section 3.4 we 
see that only the topology of the (N, L) plot really matters 
for the dependnce on coefficients Cij . 

The second simplification we make is to set v = 0. As 
mentioned, the Nodal and Lefty fronts are traveling at equal 
speed in the steady state, and are thus stationary in the co- 
moving frame. In fact, the only difference between the v — 
and u ^ cases is the steepness of their exponential profiles 
at the front. This does not qualitatively affect the results we 
present below. 

With Cnl = and v = 0, the Nodal concentration takes 
the form given by (3.8), or in dimensionless form, 



N(y) = 




and our 7-dimensional parameter space reduces to a 3-dimensional 
one, with parameters II, Cll, and L* (which is more con- 
venient for our purposes than Cln)- Our goal is to classify 
all the possible types of behavior of the Lefty profile as these 
parameters are varied. Note that Cnn and Cln are always 
greater than 1, or else N and L will never have their produc- 
tion turned on. Also, let the position where Lefty production 
turns on be y = y^', it may be on either side of the indepen- 
dent Nodal front y = 0. 

4.2 "Pinned" intervals 

It soon became apparent that the behavior of L did not only 
consist of intervals of full production sl = sql, where Cl(v) > 
1, and zero production, where Cl(v) < 1, but also intervals 
where the threshold function is pinned at CL{y) = L The 
clearest example of this is to consider what happens when 
L* < 1 as y — > — oo, where our background Nodal pro- 
file goes to a limit of N 1. Since we want L to be pro- 
duced somewhere, Cln = ClnNq > 1 necessarily (this is 
(2.6b) with L = 0). As y -> -co and N ->■ 1, L(y) = is 
mathematically inconsistent, and so Lefty production must 
turn on. Now, by definition L* is the amount of L that turns 
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off its own production when N = N n , so < 1 means 
that before L — > 1, Lefty will turn off its own production 
again. This is a "paradox" in which Lefty production cannot 
be fully on or fully off, and it is only resolved if we let the 
threshold function be pinned at 1, i.e. Cl(u — > — oo) = 1, 
such that < sz, < 1. Constraining Cl{u) means that L(y) 
is completely determined by N(y) in this pinned interval: 



i 



C LN N(y) - 1 



(4.2) 



Now, we explicitly work out the possible behaviors of L 
when L*, II and Cll are varied. L will be described in 
terms of the types of production occurring along the line, 
from — oo to oo. Full production is denoted "1", zero pro- 
duction denoted "0", and the pinned interval denoted "p". 
For example, the background Nodal profile is characterized 
as {1,0}, which is shorthand for an interval of full produc- 
tion adjoining one of zero production. In the following, we 
split up the cases into < 1 and L* > 1. We find that the 
most important parameter is the ratio of length scales, II, 
which determines how various pinned intervals arise in the 
behavior of Lefty. Note that the value of vl does not enter at 
all. 



4.3 Classification of ID model 

4.3.7 Case L* < 1: {p,0} and {p,l,0} 

As already shown, there is a pinned interval extending to 
y — > — oo, where L is completely determined. The produc- 
tion SL(y) is also completely determined, and is derived as 
follows. Substituting (4.2) back into (3.2b): 



D L d 2 N(y) 1 - „ 

-L + s L {y) = 



Cll dyi 



TL 



Now using (3.2a) to eliminate d 2 N(y)/dy 2 , Lefty produc- 
tion in the pinned interval is given by: 



C 



-(Cln[(1 - ll)N{y) + P L S N (y)] - 1) (4.3) 



LL 



where Sj^(y) 



,y>0 
,1 ,2/<0' 

The point of finding this expression is that as y —> — oo, 
we see that sl — > < 1. This simply means that as L* 
goes to zero, so does the asymptotic production of Lefty. 

The values of II, Cln or Cll in (4.3) are free to vary so 
long as L in (4.2) is between and 1 . On the other hand, the 
location of the front separating the pinned interval and the 
region of zero production must be determined from match- 
ing both L and dL j dy on either side of the boundary. It is 
possible to end up with a value of yL where Sl lies outside 
[0,1]. Indeed, checking when when this occurs using (4.3) 
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FIGURE 3 Form of wavefront when L* < 1 (section 4.3.1). 
Top: {p,0}. Bottom: {p,1,0}. Graphs are C L (y) (red dashed line), 
s L (y) (blue dot-dash line) and L(y) (black solid line); the thresh- 
old is unity. Parameter values: L* = 0.8, l L = 1 (top), l L = 1.5 
(bottom) 



demonstrates that {p,0} is inconsistent when 

LJ L > 1 



(4.4) 



Here, there are three regions instead: pinned, full produc- 
tion, and no production (in short, {p,l,0}). As LJl increases 
beyond 1, the length of the fully producing interval increases. 
Simulations and solving for the boundary conditions validate 
this (see figure 3). 

4.3.2 Case £* > 1: {1,0} and {l,p,0} 

The traveling wavefront given by (3.9) is still not guaran- 
teed when L« > 1. We can see this by considering the limit 
at which II -C 1: at the front of Lefty production, L goes 
from 1 to with a background of nearly constant N. Since 
L inhibits its own production, the apparent paradox is that 
dC 'l(u) I 'dy > 0, i.e. the regions of zero and full production 
are reversed! Simulations show that as II is decreased and 
the slope of L made steeper, once the slope exceeds that of 
the "pinned" Lefty (4.2), {1,0} will transition to {l,p,0}. Ba- 
sically, requiring that dCL (y) /dy < at the front is equiva- 
lent to this condition. 

To write out explicit conditions for entering the {l,p,0} 
state, we should consider both cases yL > and yL < 0, 
for which N takes different forms given in (4.1). Thus, the 
condition for {l,p,0} is that Cx,(y) crosses the threshold in 



Biophysical Journal: Biophysical Letters 



L08 



Biophysical Journal: Biophysical Letters 



0.5 




0.5 




FIGURE 4 Form of wavef ront when L* > 1 (section 4.3.2). Top: 
{1,0}. Bottom: {1,p,0}. The curves C L (y) (red dashed), s L (y) 
(blue dot-dash) and L(y) (black solid) as in Figure 3. Parameter 
values: L» = 1.5, l L = 1 (top), l L = 0.4 (bottom) 



the other direction: dCx(yi)/dy > 0, giving: 

2 



y L > : ll 1 > 1 . 
y L < : IT 1 > 2L* - 1. 



(4.5a) 
(4.5b) 



(Note that both conditions require II < 1, i.e. the Nodal 
length scale is greater than Lefty, as argued above in the lim- 
iting case.) If the inequalities are in the opposite direction, 
then we are back to the traveling wavefront given by (3.9), 
or in our notation 1,0. Here we can solve for the Lefty con- 
centration profile fully: using (3.10) we find that Ay > 
when 

\c LN - \c LL < 1. (4.6) 



5. Solutions oscillating in time 

In the previous section(s), we explored the limiting case of an 
autonomous, Lefty-independent Nodal production. We now 
address new phenomena which are possible when Nodal is 
significantly affected by Lefty. The most striking of these is 
the possibility of a time-oscillating solution, which we will 
study in the special case that all concentrations are uniform 
in space. The experimentally pertinent motivation is whether 



Interval patterns 


Parameter conditions 


{1,0} 


L* > 1, dC L (y L )/dy < 


{1,P,0} 


L. > 1, dC L (y L )/dy > 


{p,o} 


<C 1, > 


{p,l,0} 


<C I, <C 



200 

y 



200 

y 



Table 1 Conditions for possible Lefty interval patterns, given 
a Nodal wavefront for Cnl = and v = 0. Note that there are 3 
independent parameters: l L , C LL , and L» (or C LN ) 



it is possible to generate traveling solutions, consisting of 
pulses of Nodal and Lefty, periodic in space and in time; 
one expects such solutions to be possible only if uniform 
oscillations are possible in the model. 

By examining the geometry of the (N, L) plane [Figure ??], 
it can be seen that only case (d) has the possibility of os- 
cillations; this is shown with more detail in Fig. 5(b). Note 
that the fixed point in the center, since it has Np < Nq and 
hp < Lq, corresponds to production rates less than satu- 
ration (sjy < snq and sl < slo), so we are in a "pinned" 
regime for both N and L, which is self-consistent with being 
on the isoclines for both. 



5.1 Parts of each cycle 

The period of an oscillation must consist of four phases [as 
shown in Fig. 5(a)]: 

(I) Nodal and Lefty both on; Lefty grows faster than Nodal, 

until Nodal turns off; 
(II) Nodal off and Lefty on; Nodal decays, until Lefty turns 
off; 

(III) Nodal and Lefty both off; Lefty decays faster than Nodal, 
and they remaining Nodal that has not yet decayed is suf- 
ficient to turn on Nodal again; 

(IV) Nodal on and Lefty off; Nodal grows, until it is sufficient 
to turn Lefty on too (and back to phase 1). 

Evidently, a key condition is on the degradation times, which 
govern the growth rates of Nodal or Lefty : we need 



TL = T L /T N < 1. 



(5.1) 



(Of course, the actual magnitude of the degradation times 
simply sets an overall time scale.) 

We can describe the history with a series of times ti, i=l,2 

at each of which one of the signals turns on or off. Between 
these times, the production rates are constant (either zero or 
saturated), and the differential equation is solved by a simple 
exponential: 



N(t) = 



N{U)e 



-{t-n)/T N 



1 - [1 - N(U)]e-^- u ^ TN 
and similarly for L(t). 



sjv = 
sn 



snq, 
(5.2) 
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To visualize the dynamics, it is convenient to draw the 
trajectory of (N(t), L(t)) in the (N, L) plane; in our equa- 
tions, the time derivatives N = dN/dt and L = dL/dt 
are fuctions only of (N,L). The lines Cn(N,L) = and 
Cl(N,L) — on this plot are the isoclines, meaning the 

places where (respectively) N and L change sign. (Isoclines 
were also used in the analysis of Ref. [13].) Each time the 
trajectory intersects an isocline is one of the times U; let the 
concentrations at these times be Ni = N(U), Li = L(U). 
We can find the trajectory curve, and thus find solutions of 
the dynamical equations, by eliminating time and consider- 
ing only the discrete map (Ni, Li) — > (Ni+i , Li+\ ). We next 
find the functional form of this map. 

Consider phase III of the cycle, during which [according 
to (5.2)] N(t) = A^e-t*"**)/™ and L{t) = L i e-^- t ^l TL . 
Eliminating time, we get 



(a) 



L 



- 1 

N 



(5.3a) 



The point (iVj+i , Li+i) is the intersection of the curve (5.3a) 

0. In phase IV, the trajectory 



with the isocline Cn(N, L) 
would be 



L 



1-N 
l-Ni 



(5.3b) 



and similarly in the other four phases of the cycle. Notice 
that the curve can possibly intersect the isocline in phases 
I or III only if fx < 1. (For example, whehn fx = 1 the 
trajectory in each phase is a straight line to one corner of the 
square.) 

To give the cycling behavior, it is necessary (but not suf- 
ficient) that the isoclines intersect as shown in Figure 2(b), 

which depends on two inequalities. The N and L isoclines' 
respective intercepts at L = satisfy 1/Cnn < 1/Cxjv; 
their intercepts at L = 1 must be in the reverse order, (1 + 
Cnl) I Cnn > (1 + Cll) /Cln- The two inequalities can 
be written together as 



Cll + 1 Cln j 



Cnl + 1 Cnn 



(5.4) 



This tends to be satisfied When conditions (5.4) are satisfied, 
the isoclines cross at a fixed point (iV*, L„) where 



N* = — 



Cnl — Cll 



L, 



ClnCnl — CnnCll 

Cnn — Cln 
ClnCnl — CnnCll 



(5.5a) 
(5.5b) 



5.2 Stability conditions 

The trajectory always tends to spiral around (iV«, L t ). as il- 
lustrated in Figure 5(b). However, there are different pos- 
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FIGURE 5 (a). Dynamics of the concentrations N(t) and L(t) 
of Nodal and Lefty, and their production rates s N (t) and s L (t), in 
the case of periodic oscillation (schematic). Phases of the cycle 
are labeled I, II, III, and IV. (b). Phase plane of the concentrations. 
Heavy lines indicate the sign change of N(N, L) (green) and of 
L(N, L) (red). A trajectory (blue) is shown, which (initially) is 
spiraling out from an unstable fixed point (N„,L t ). Filled and 
open circles are stable and unstable fixed points. The trajectory 
from the heavy cross is tangent to the next isocline; any trajec- 
tory starting higher goes to (0, 0). 

sibile limiting behaviors depending on whether the trajec- 
tory spirals inwards (stable fixed point) or outwards, and on 
where it ends up. There is always a stable fixed point at 
(N,L) = (0,0) [the empty system] and provided the iso- 
clines intercept the upper edge of (N, L) space [as in Figure 
5(b)], there is also a stable fixed point at (N,L) = (1, 1) 
[both Nodal and Lefty turned on and saturated]. Some of the 
possibilities: 

(1) The fixed point is stable; there is one unstable cycle, such 
that a trajectory starting inside it tends to (iV*, L*) and a 
trajectory starting outside it tends to (0, 0) or (1, 1). 
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(2) There are no cycles (stable or unstable); the fixed point is 
unstable, and spirals out to (0, 0) or (1, 1). 

(3) The fixed point is unstable, and spirals out to a stable cy- 
cle (beyond which is an untable cycle, as in case (1)); this 
is the case of interest to us. 

To evaluate the stability, we must linearize near the fixed 
point in terms of SN = N — N*,6L = L—L*. The trajectory 
(e.g. Eqs. (5.3a) or (5.3b)) becomes 8L = m p h5N, also 
rriph is the slope in each respective phase of the cycle, ph — > 
I,II,HI, or IV. Thus, 



1 L* 
mm = — -~- 



1 / L t 



mrv = 



(5.6a) 



(5.6b) 



tl V 1 - N* 

and similarly for the rest. A convenient viewpoint on the con- 
dition is that 

"tph = f L m phfi (5.7) 

where m p h,o is the slope of the line from (iV*,Z*) to the 
appropriate corner of the square. 

Meanwhile, the slopes of the N = and L = isoclines 



mjy = Cnn/Cnl, 
m L = C LN /C LL . 



(5.8a) 
(5.8b) 



Evidently, the additional (and sufficient) condition to get a 
spiraling behavior is that mjv < mm, mi < mL- Since 
tl < 1 and m p h is of order 1/t£ for ph = I, II, III, IV, 
it follows that is typically large and hence Cll must be 
relatively small. 

Solving one linear equation for each phase of the cycle, we 
find after four isocline intersections that (8N i+ 4, 8L i+ 4) = 
A(5Ni, SLi), where 



A = 



i?i R 



in 



RjiRtv 



where 



R 



P h 



(5.9) 



(5.10) 



niN - m p h 

for ph = I, II, III, IV. Thus, the fixed point (iV*,Z*) is sta- 
ble if and only if A < 1. An interesting special case is when 
(iV*,Z*) = (1/2, 1/2): in that case, A = 1 so the fixed point 
is always marginal. 

To get a better understanding of the stability, consider the 
case that mjv <C mi, miv <C mi. Then 



1? 



ph 



1 



m ph \ m L mph 
Noting that mimin/mnmiv = 1, we get 
1 



(5.11) 



Awl 



TLm L 



mi , + mm n - mn,o - mrv.o + 



+T L m N I jtij q 1 + m n j - TOjj q - mjy ) . (5.12) 



The sums in parentheses do not have a factor of fx; they de- 
pend only on (TV*, L„). Both sums tend to be positive when 
the fixed point is farther from the line L = N than it is from 
the line L = 1 — N, or both negative in the opposite situa- 
tion. Thus, in the former situation, the fixed point tends to be 
unstable when 

(m L m N )^ 2 < f L . (5.13) 

In the opposite situation, the fixed point tends to be unstable 
when the inequality is in the other direction. 

The oscillation period T can be inferred from the orbit on 
the (N, L) plane by going back to equation (5.2). Note that 

N and L do not go to zero approaching the isocline line, but 
rather they approach a positive or negative constant coming 
from the respective sides. In other words, the four legs of 
the cycle on the (TV, L) plane are each traversed at a roughly 
constant speed. Hence a cycle that forms a small loop around 
(iV*,-L*), as happens just beyond the instability, has its pe- 
riod T proportional to its oscillation amplitude. 

It follows that if the parameters are close to the instability 
threshold, we can construct a spatially periodic solution of 
intervals containing oscillations which are shifted in phase 
with respect to each other, provided that the spatial period is 
large compared to the distance that either signal can diffuse 
in time T: locally, this situation is equivalent to the uniform 
one, since the signal from contrasting intervals does not have 
time to propagate. It follows that periodic traveling waves 
are also possible, as well as a standing wave made up of al- 
ternating domains that each oscillate similar to the uniform 
oscillations, but in opposite phase to each other. 

6. Conclusion 

In this work, we first sketched (Sec. 2) the realistic two- 
dimensional model developed in [14], which exhibits a uni- 
formly moving posterior-to-anterior front of Nodal produc- 
tion on one of the two stripes representing the lateral plate 
mesoderm, followed by a front of Lefty production on the 
midline, just like the experimental observation [11], and which 
is known to propagate the left/right asymmetry from the tail- 
buds through the entire embryo. We then set up a version 
of this model in one dimension (Sec. 3), representing the 
anterior-to-posterior axis, which captures most properties of 
front propagration. The one-dimensional space represents a 
combination of the left LPM and the midline, and its main 
qualitative defect compared to the two-dimensional system 
is that it misses the lag time for the Nodal signal N to dif- 
fuse from the LPM to the midline, or for the Lefty signal 
L to diffuse from the midline to the LPM. The key param- 
eters of the model were (1) coefficients [Sec. 2.2] quanti- 
fying how much Nodal promotes the two signals, and how 
much Lefty inhibits them, and (2) different diffusion con- 
stants and/or degradation times for the two signals, allowing 
us to define a typical length scale for either one, over which 
its concentration varies in front of and behind a sharp spatial 



Biophysical Journal: Biophysical Letters 



L011 



Biophysical Journal: Biophysical Letters 



step of the production; all in all, we found seven dimension- 
less ratios, each of which is a nontrivial parameter (Sec. 3.3). 

We then set out to classify, analytically, the possible be- 
haviors of the model in a steady state, and see how the de- 
pended on the parameter values. A key aid in identifying the 
regime is make a two dimensional plot (Sec. 3.4) of how 
the homogeneous concentrations (N,L) would evolve: the 
very existence of a front depends on having a fixed point on 
this plane representing nonzero Nodal production, alongside 
the fixed point of zero production (which is always present). 
Our main focus was on uniformly moving fronts, but in fact 
one can specialize (without real loss of generality, and with 
considerable mathematical simplification) to the special case 
of zero front velocity [Sec. 3.2]. The key special feature we 
found was the possibility that, instead of showing a sharp 
front where the production rapidly switches from zero to 
maximum, the Lefty production (on the midline in the real 
embryo) might be varying continuously over an extended in- 
terval (Sec. 4). This was due to a sort of feedback or buffer- 
ing, in which the Nodal and Lefty concentrations compen- 
sate exactly such that every point along this interval is right 
at the threshold for Lefty production ("pinned"). 

Finally, in Sec. 5 we changed focus to the possibility of 
a temporally periodic behavior; this parameter regime al- 
lows both the possibility of a steady limiting state which is 
"pinned" (at threshold) for both Nodal and Lefty. One may 
get oscillations if Lefty degrades much faster than Nodal, 
and we characterized the exact mathematical criterion for the 
steady state to be unstable to developing oscillations. 

Having found the relationship of these behaviors to the pa- 
rameter values, in future work we will be able to find the cor- 
responding parameter regimes in the realistic two-dimensional 
model, with the ultimate goal of finding a set that gives the 
observed phenotypes. We obtained many behaviors that are 
not seen experimentally. (It is still uncertain whether the Lefty 
production has a sharp front or an extended one like the 
"pinned" case.) This is a positive feature of the model, for 
it will allow the exclusion of whole domains of parameter 
values. 

Our results might motivate further studies of the possible 
behaviors of propagating L/R asymmetry in the zebrafish and 
other embryos. The pertinent measurements in real systems 
that relate to phenomena in this model are 

(1) The measured front velocity constrains the velocity scale 
v N . 

(2) The sharpness of the N and L concentration fronts gives 
the length scales In and II and/or indicates the possibility 
of the "pinned" interval; however, experiments that mea- 
sure mRNA are in effect probing the production rate and 
do not access In and I l directly. 

(3) The lag between the Nodal front and the Lefty front, gives 
information on the difference between Nodal and Lefty 
threshold parameters. Experimentally the Lefty front lags, 
but many parameter sets give a Lefty front that spatially 



leads (even though causally it is downstream from Nodal), 
so observation can limit parameters to a certain window. 

In addition, they suggest how changing a parameter (due to a 
mutation or some external perturbation) could give an exotic 
result, such as non-monotonic front. 

What the present study did not address is how the Nodal 
production is first initiated, and why it only happens along 
the left LPM. That depends on the two-dimensional geome- 
try, e.g. the relative separations of the tailbud from the foot 
of the LPM stripes and the midline stripe of Oep-producing 
cells. Furthermore, it is possible the real system could propa- 
gate, e.g., a finite pulse of Nodal that terminates, but that will 
be seen only if such a signal is produced at the posterior end. 
These questions must be left to a proper two-dimensional 
study. 
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